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Analytic Properties and Covariance Functions 
of a New Class of Generalized Gibbs Random 

Fields 

Dionissios T. Hristopulos and Samuel N. Elogne 
Abstract 

Spartan Spatial Random Fields (SSRFs) are generalized Gibbs random fields, equipped with a coarse-graining 
kernel that acts as a low-pass filter for the fluctuations. SSRFs are defined by means of physically motivated spatial 
interactions and a small set of free parameters (interaction couplings). This paper focuses on the FGC-SSRF model, 
which is defined on the Euclidean space R'' by means of interactions proportional to the squares of the field 
realizations, as well as their gradient and curvature. The permissibility criteria of FGC-SSRFs are extended by 
considering the impact of a finite-bandwidth kernel. It is proved that the FGC-SSRFs are almost surely differentiable 
in the case of finite bandwidth. Asymptotic explicit expressions for the Spartan covariance function are derived for 
d = 1 and d — 3; both known and new covariance functions are obtained depending on the value of the FGC-SSRF 
shape parameter. Nonlinear dependence of the covariance integral scale on the FGC-SSRF characteristic length is 
established, and it is shown that the relation becomes linear asymptotically. The results presented in this paper are 
useful in random field parameter inference, as well as in spatial interpolation of irregularly-spaced samples. 

Index Terms 

parameter inference, Geostatistics, Gaussian, correlations. 

I. Introduction 

Spatial Random Fields (SRF's) have a wide range of applications in hydrological models [14], [26], [30], 
petroleum engineering [16], environmental data analysis [10], [25], [36], mining exploration and mineral reserves 
estimation [4], [15], environmental health [11], image analysis [40], [41], medical image registration [3] and brain 
research [9], [28], [34] among other fields. 

A spatial random field {X{s,uj) e M; s e D{L) C M'^;a; G il} is defined as a mapping from the probability 
space (il, A, P) into the space of real numbers so that for each fixed s, X{s, lo) is a measurable function of u) [1, p. 
3]. D{L) is the domain within which the SRF is defined and L is a characteristic domain length. An SRF involves 
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by definition many possible states [10, p. 27], [42], denoted by u). In the following, for notational simplicity we 
suppress the dependence on u). 

Realization of a particular state is determined from a joint probabihty density function (p.d.f.) /x [-''^(s)]. The 
p.d.f. depends on the spatial configuration of the field's point values. For spatial data the term sample refers to N 
values X{si) from a particular state at the measurement locations {sj, i = 1, . . . , N}, representing a single state 
of the SRF (an observed realization). 

For irregularly-spaced samples, practical applications involve determining the statistical parameters of spatial 
dependence and interpolating the data on a regular grid. An observed realization, X*(s), can be decomposed into 
a deterministic trend m^{s), a correlated fluctuation SRF Xx(s), and a random noise term, e(s), i.e., X*{s) = 
mx(s) + X;^(s) + e(s). The trend is a non- stationary component representing large-scale, deterministic variations, 
which presumably correspond to the ensemble average of the SRF, i.e. m^{s) — E[X{s)]. 

The fluctuation corresponds to variations that involve smaller spatial scales than the trend. It is assumed that 
resolvable fluctuations exceed the spatial resolution A. The component e(s) represents inherent variability below 
the resolution cutoff, or completely random variability due to uncorrelated measurement errors. It will be assumed 
that e(s) is statistically independent of the SRF Xx{s), and can be treated as Gaussian white noise. 

It is assumed in the following that the trend has been estimated and removed. This work focuses on modeling 
the correlated fluctuations. 

For many geostatistical apphcations, the fluctuation can be viewed as a weakly stationary SRF, or an intrinsic 
random field with second-order stationary increments [42, pp. 308-438], [29]. An SRF is weakly stationary if its 
expectation m^{s) is independent of the location, and its covariance function Gx(s,s -|- r) depends only on the 
spatial lag r. In the following, the term 'stationary' wiU refer to weak stationarity. 

Furthermore, a stationary SRF is statistically isotropic if the covariance function depends only on the Euclidean 
distance between points but not on the direction of the lag vector, i.e., Gx(|r|), where |r|) is the Euclidean norm 
of the vector r. The isotropic assumption is not restrictive, since the anisotropic parameters can be inferred from 
the data and isotropy can be restored by rotation and rescaling transformations [13], [17], [20], [21]. 

For isotropic, short-ranged SRF's, one can define a single integral scale [30, p. 22], given by the integral of 
the covariance function along any direction in space. Since the parameters of the fluctuation SRF are determined 
from the available sample by employing the ergodic hypothesis [42, p. 29], [27, p. 30], the integral scale must be 
considerably smaller than the domain size L. 

The distribution of Gibbs random fields is expressed in terms of an energy functional H[X\{s); 0], where 6 is 
a set of model parameters as foUows: [41, p. 51] 



The constant Z{0), called the partition function normahzes the p.d.f. and is obtained by integrating exp {—H[Xx{s) ; 6]} 
over aU the reaUzations. 




(1) 
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Gaussian SRF's used in classical geostatistics can be included in the formalism of Gibbs SRF's if the energy 
functional is expressed as follows: 



where the precision matrix, [G] ~ {6) is the inverse of the covariance matrix; the latter is determined directly from 
the data by fitting to parametric models. Note that in Eq. (|2jl there is no explicit resolution scale. 

The rest of this paper is structured as follows: Section (|n|i gives a brief overview of Spartan spatial random fields. 
Section JIIH focuses on general (for any d) properties of the covariance function. Section iTVi proves the property of 
sample differentiability for a specific class of SSRF models. Section (|3 focuses on a one-dimensional SSRF model 
and obtains explicit expressions for the variance and the integral scale, as well as expressions for the covariance 
function valid in the infinite-band limit. Section dVH obtains explicit respective expressions in three dimensions. 
Section dVIH summarizes the contributions derived in this paper Finally, the calculations used in deriving the 
variance (finite-band case) and the covariance functions (infinite-band limit) are presented in detail in Appendices. 



The term Spartan indicates parametrically compact models that involve a small number of parameters. The Gibbs 
property stems from the fact that the joint probability density function (p.d.f.) is expressed in terms of an energy 
functional H[X{s)], i.e., /x[X(s)] = exp{—H[X{s)]}. Use of an energy functional containing terms with a clear 
physical interpretation permits inference of the model parameters based on matching respective sample constraints 
with their ensemble values [18]. Thus, the spatial continuity properties can be determined without recourse to the 
variogram function. 

Estimation of the experimental variogram, using the classical method of moments or the robust estimator, [12], 
involves various empirical assumptions (such as choice of lag classes, minimum number of pairs per class, lag and 
angle tolerance, etc. [15, pp. 75-123], [38, pp. 44-65]). In addition, inference of the 'optimal' theoretical model 
from the experimental variogram presents considerable uncertainties. For example, least-squares fitting may lead to 
sub-optimal models. This is due to the proportionally larger influence of larger lag distances that correspond to less 
correlated fluctuations. This situation often forces practitioners to search for a visually optimal fit [38, pp. 48-49] 
of the experimental variogram with a model that yields better agreement in the short distance regime. It may be 
possible to address some of these shortcomings more effectively in the SSRF framework. 

The recently proposed FGC-SSRF models [18] belong in the class of Gaussian Gibbs Markov random fields. The 
Markov property stems from the short range of the interactions in the energy functional. However, the development 
of SSRFs does not foUow the general formalism of GMRF's [5], [12], [31], [32], [41]. The spatial structure of die 
SSRFs is determined from the energy functional, instead of using a transition matrix, or the conditional probability 
(at one site given the values of its neighbors). Model parameter inference focuses on determining the coupling 
strengths in the energy functional, instead of the transition matrix. The GMRF formalism and related Markov 
chain Monte Carlo methods can prove useful in conditional simulations of SSRFs. Methods for the non-constrained 




(2) 
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simulation of SSRFs on square grids and two-dimensional irregular meshes are presented in [19], [22]. In fact, 
SSRF models without the Gaussian or Markov properties can be constructed. 

The energy functional of SSRFs involves derivatives (in the continuum), suitably defined differences or kernel 
averages (on regular lattices and irregular networks) of the sample [35]. In all cases, the energy terms correspond to 
identifiable properties (i.e., gradients, curvature). As we show below, in the continuum case the energy functional 
is properly defined only if X\{s) involves an intrinsic resolution parameter 'A'. This parameter is introduced by 
means of a coarse-graining kernel. Hence, formally, the SSRFs belong in the class of generalized random fields [2, 
p. 44], [42, pp. 431-446]. 

The resolution limit A is a meaningful parameter, since a spatial model can not be expected to hold at every 
length scale; also, in practice very small length scales can not be probed. In contrast, classical SRF representations 
do not incorporate a similar parameter. The resolution limit imphes that the covariance spectral density acts as a 
low-pass filter for the fluctuations. 

For lattice-based numerical simulations, the lattice spacing provides a lower bound for A. The latter, in frequency 
space corresponds to the upper edge, l/2a, of the Nyquist band; equivalently, in wavevector space it corresponds to 
the upper edge, n/a, of the first Brillouin zone. For irregularly spaced samples there is no obvious cutoff a priori; 
hence the cutoff is treated an a model parameter to be determined from the data. Consequently, we use the spectral 
band cutoff, kc, as the SSRF model parameter instead of A. 

A. The FGC Energy Functional 

A specific type of SSRF, the fluctuation-gradient-curvature (FGC) model was introduced and studied in [18]. 
The FGC energy functional involves three energy terms that measure the square of the fluctuations, as well as their 
gradient and curvature. 

The p.d.f. of the FGC model involves the parameter set 6 — {ijo, ^, kc): the scale coefficient 770 determines the 
variance, the shape coefficient 771, determines the shape of the covariance function, and the characteristic length ^ 
is Unked to the range of spatial dependence; the wavevector kc determines the bandwidth of the covariance spectral 
density. If the latter is band-limited, kc represents the band cutoff and is related to the resolution scale by means 
of kcA w 1. 

More precisely, the FGC p.d.f. in M'' is determined from the equation: 



where 6' = (771, ^,kc) is the reduced parameter set, and iifgc is the normalized (to 770 = 1) local energy at s. The 
functional hfgc [X;^(s); 0'] is given by the following expression 



The explicit, non-linear dependence of Eqs. (|3} and @ on the parameters 77o,77i,^,kc that have an identifiable 




(3) 



/7fgc [Xx{sy,e'] 



(4) 
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physical meaning is preferable to the use of linear coefficients for the variance, gradient and curvature, because it 
simplifies the parameter inference problem and allows intuitive initial guesses for the parameters. 

III. The FGC Covariance Function 

The FGC energy functional has a particularly simple expression in Fourier space. Let the Fourier transform of 
the covariance function in wavevector space be defined by means of 

G,;A(k;0) = y"dre-^''"-Gx;A(r;0), (5) 

and the inverse Fourier transform by means of the integral 

G,;A(r;0)-^^ ydke^'^ "-G,;A(k;0). (6) 

In the following, we suppress the dependence on 9 when economy of space requires it. 
The energy functional in Fourier space is given by 

i7fge[XA(s);0] - A_.X^{\^) [G.a]-'(k)XA(-k). (7) 

Note that the interaction is diagonal in Fourier space, i.e., the precision matrix [Gx;A]^(k) couples only components 
with the same wavevector value. 

Also, for a real-valued SSRF X{s) it follows that Xa(— k) = x\{k). For a non-negative covariance spectral 
density, it follows from Q that the energy is also a non-negative functional. 

The covariance spectral density follows from the explicit expression: 



Gx;A(k) = /,(k;0 ) QA(k) 



2 



(8) 



where = (?7o,?7i,C)' Q\(^) is the Fourier transform of the coarse-graining kernel and 

In [18], [19], a kernel with an isotropic boxcar spectral density, i.e, with a sharp wavevector cut-off at kc, was 
used. The boxcar kernel leads to a band-limited covariance spectral density Gx;A(k). This kernel will be used here 
as well. It involves a single parameter, i.e., kc, which facilitates the inference process. Nonetheless, it is not the 
only possibility. 

For this functional to be permissible, the covariance function must be positive definite. If kc^ is considered as 
practically infinite, application of Bochner's theorem [6], [42, p. 106], permissibility requires r/i > —2, as shown 
in [18]. For negative values of 771 the spectral density develops a sharp peak. Gx;A(k) tends to become singular 
as rji approaches the permissibility bound of —2. In early investigations [18], [37], kc was treated as an a priori 
known constant so that kc^ >> 1. However, it is also possible to infer the value of kc from the data [35]. In this 
case, the permissibility criterion is modified as follows; 
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Theorem 1 (Permissibility of FGC-SSRF): The FGC-SSRF is permissible (i) for any kc if r/i > —2 and (ii) for 
7]! < -2, provided that kc^ < -^y/lvil - where A = jy^i^ - 4|^. 

Proof: Let us assume that = {fc e M : Qxik) 0}. Based on Eq. we obtain /^(fc; S") = 770 ^'^/n(fc^), 
where n(x) = 1 + r/ix^ + x"*. Then, n(x) = (x^ — yi){x^ — 1/2), where j/1.2 = {—ill i ^)/2- Bochner's theorem 
requires that n(fc^) > 0, Vfc G 1?^. The case for 771 > —2 is proved in [18]. For r/i < —2 it follows that 1/12 — 



{\r]i \ ± A)/2 > 0. Hence, Bochner's theorem is satisfied if Vfc e Dfc : fc^ < min(^/yr, ^/yj) = a/(|77i| - A)/2. ■ 
Remark 1: Bochner's theorem is also satisfied if Vfc e Dk : fc^ > max(y^, = y^|^i|+A)/2. This case 
corresponds to a coarse-graining kernel that acts as a high-pass filter, and is not relevant for our purposes. 

The spectral representation of the covariance function is given by means of the following one-dimensional 
integral, where Jd/2-i{i') is the Bessel function of the first kind of order d/2 — 1, 

G^-Mr}-j^^J^ ^S + ^^(fc^)2 + (fc^)4- (10) 
In Eq. (llOt and in the following, we take r and k to represent respectively the Euclidean norms of the vector r 
and k. Only in d = 1, we will use \r\ and |fc| to denote the norm (absolute value). The Bessel function can be 
expanded in a series as follows [39, p. 359], where T{x) is the Gamma function: 

A. The Variance 

We investigate the dependence of the variance, = Gx;a(0) on the SSRF parameters. The covariance function is 
well behaved at zero distance, in spite of the singular factor r''/^^^ dividing the integral in Eq. ilO\ . This singularity 
is canceled by the leading-order term of Jd/2~i{kr) as r — > 0, which is given by Jd/2-i{kr) ^ (1^)'^^^^^ /T{d/2). 
For r = only the leading-order term of the expansion il 1> contributes. Hence, we obtain 

kc 



a. 



2 



770^'^ f dk k 



2(d-2)/2r(|) (2^)rf/2 J l+^i(fc^)2 + (fc^)4' 



and using the variable transformation x — k£,, it follows that: 

2(<i-2)/2 r (I) (27r)'i/2 y l + 7yix2+x4- ^ ' 

This integral can be explicitly evaluated for any d as a function of 771 and kc In the infinite-band case (kc £, 00), 
the variance integral exists only for d < 4. 

B. The Integral Scale 

The integral scale of the covariance function for isotropic SRFs is given by the equation: 



G.;a(0) 



(13) 
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Using Eq. (|5} with k = 0, and Eq. (|9} for the DC component of the spectral density, the integral scale follows from 



W) = 



G.a(0) 



'70 



G.a(0) 



(14) 



In sections (jVj and (IVI> we derive explicit expressions for the variance and the integral scale in d = 1 and d = 3, 
and we study their dependence on 771 and kc ^. These expressions show that the integral scale in the preasymptotic 
regime is a nonlinear function of the characteristic length, in contrast with most classical covariance models. 

IV. Existence of FGC-SSRF Derivatives 

In this section we prove that the band-limited FGC-SSRF models have differentiable sample paths with probabiUty 
one. Conversely, for the infinite-band case only the first derivative of the SRF exists in li = 1. 

Many of the covariance models used in geostatistics are non-differentiable (e.g., the exponential, spherical, and 
logistic models). Notable exceptions are the Gaussian model (which leads to very smooth SRF reaUzations) and the 
Whittle-Matern class of covariance functions [12], [33]; the latter include a parameter that adjusts the smoothness 
of the SRF. Hence, band-limited SSRFs enlarge the class of available differentiable SRF models. 

Non-differentiable covariance models are often selected for processes the dynamical equations of which are not 
fully known or can not be solved, based solely on the goodness of their fit to the experimental variogram. However, 
this does not imply that the sampled process is inherently non-differentiable. If most of the candidate models are 
non-differentiable, or very smooth differentiable ones, it is not surprising that the former perform better than the 
latter. Differentiable SRF models with controlled roughness may provide equally good candidates. 

A. Partial Derivatives in the Mean Square Sense 

The existence of first and second order derivatives of X\ (s) in the mean square sense is necessary to properly 
define the FGC-SSRF. This follows since the energy functional, given by Eq. @, involves the spatial integral of the 
squares of the gradient and the Laplacian. Assuming ergodicity, these integrals can be replaced by the respective 
ensemble mean multiplied by the domain volume (in W^). 

For stationary Gaussian SRFs, a sufficient condition for the field partial derivatives to exist in the mean square 
sense [1, p. 24] is the following: 

Let it — {nil , . . . , Ud) be a vector of integer values, such that ni + . . . + nd = n. The nth-order partial derivative 
9"X(s)/(9s"^ . . . ds^^'' exists in the mean square sense if the following derivative of the covariance function exists 
[2] 



G(")(0) 



2ni 



. dr 



(15) 



r=0 



Theorem 2 {Mean-Square Differentiability): For FGC Spartan Spatial Random Fields with a band-limited covari- 
ance spectral density, the partial derivatives of any integer order n are well defined in the mean square sense. 
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Proof: The FGC SSRFs are stationary and jointly Gaussian. Hence, the existence of the covariance partial 
derivative M5\ needs to be proved. It suffices to prove that |Gi"''(0)| exists. Equivalently, it suffices to prove the 
convergence of the following Fourier integral: 



|Gi"Ho)| = m^'' j and j dk Qx{k) 



. 2ni 1,2nd ud-l 



2 



(16) 



where dVld is the solid angle differential. Let ki = k cos (pi and then define the following integral over the unit 
sphere: Zd — J dfldcoscpi . . .cos^^. Note that Zd < Sd = J dfld = 2i:'^/'^ /T{d/2), where 5^ is the surface area 



of the unit sphere in d dimensions. Then, the |Gi"''(0)| is given by: 



\G^J^\n^meZd / dk 



(17) 





If Q\{k) is the boxcar kernel, [^"^(O)! in Eq. ([^ is expressed in terms of the following integral: 

1 + ryi x2 + X- 

The integral on the right-hand side on the inequality (I18> converges for all d and n. This establishes the sufficient 



|G<"nO)l=^or^"^rf I d^ ^ _ . (18) 







condition for the existence of partial derivatives in the mean-square sense. ■ 

Remark 2: The proof focused on the boxcar kernel, but the same arguments can be used for any kernel that 
decays at large k faster than a polynomial. 

If is fixed, the integral of Eq. ( I18> is proportional to implying that the SSRF is smoother for larger 

^. For X >> 1 the integrand behaves as x''+'^"^^. If ^ is fixed, the contribution of the large x in the integral of 
Eq. ( I18> scales as (kc)^" (kc^*^ ^- This scaling implies that the roughness of the SSRF increases with kc. 

Corollary 1 (Infinite-Band Case): For FGC SSRFs with an infinite band, only the second-order partial derivative 
of the covariance exists in c? = 1. Higher-order derivatives do not exist even in d = 1, and derivatives of any order 
are not permissible in any d> 1. 

Proof: For the boxcar kernel with kc — > oo (i.e. in the absence of smoothing), the integral in the right-hand 
side of Eq. ( I18t converges for d + 2?! < 4 and diverges in all other cases. Convergence is attained only for d < A 
and n = or for d — \ and n = 1. Hence, only the first-order derivative in d = 1 exists in the mean-square sense. 
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Remark 3: If a kernel that behaves asymptotically as 
convergence condition becomes d + 2n < 4 + p. 



Qx{k) 



IX k P is used instead of the boxcar, the 



The existence of the first-order derivative is not sufficient to guarantee that the FGC SSRF has second-order 
derivatives in the mean square sense. Hence, a band limit is necessary to obtain well defined second derivatives 
and the energy functional of Eqs. Q and (0). 

Below, we derive explicit asymptotic expressions for the co variance function in d — 1,3 that do not admit second- 
order derivatives. These should be viewed as limit forms of the FGC-SSRF model when kc oo. However, it is 
also shown that the asymptotic expressions are accurate estimators of the covariance function for any kc, provided 
that kc^ > Vd, where Vd is a dimension-dependent constant. 

B. Differentiability of Sample Paths 

The existence of differentiable sample paths presupposes the existence of the partial derivatives in the mean 
square sense. In addition, a constraint on the rate of increase of the negative covariance Hessian tensor near the 
origin must be satisfied [1, p. 25] to ensure the existence of derivatives with probability one. 

The negative covariance Hessian tensor is defined as follows: 

G.;.,(r)-- ^^^^^^ . (19) 

The sufficient condition for the existence of the derivative diX{s) requires that for any < r < 1 there exist 
positive constants and e^, such that the following inequality is satisfied: 

G,;,. (0) - G,;,, (r) < \ . (20) 

Theorem 3 (Existence of Path Derivatives): The FGC-SSRFs with a band-limited covariance spectral density 
have differentiable sample paths. 

Proof The FGC SSRFs are jointly Gaussian, stationary and isotropic random fields. For an isotropic SRF, 
the value of the partial derivative of G^{r) is independent of direction. Therefore, Gx-ii{r) = —AGx{r)/d, where 
AGx(r) = ^X]i'=iGx;M(r) is the Laplacian. Hence, it is sufficient to prove the validity of the inequality ( I20> for 
the Laplacian, i.e., 

- [AG.(O) - AG.(r)] < (21) 
Let us define the following function: 

Cx(r) - - [AG^O) - AGx(r)] . (22) 
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In light of Cx(r), the sufficient condition Mil becomes: 

Cx(r) < (23) 

For r — *■ the right hand side in the inequality \23\ tends to zero. Hence, the sufficient condition requires Cx(0) < 0. 
This is satisfied since Cx(0) = as it follows from the definition (I22> . 
For < r < 1, the condition can be expressed as 

Cx(r)|logr|i+^ < c. (24) 
The Laplacian of the covariance function is given by the following integral in wavevector space: 

{27r)d/^ r^/^-^ J l + ViikO^ + ikO^' 



The Bessel function in Eq. i25\ is expanded using the series The nth-order term in the expansion is cx 

j,^2n+d/2-i ^ jjjug canceling the r'^/^^'^ dependence in the denominator of — AGx(r). The leading (?i ~ 0) term 
of — AGx(r) is independent of r, while all other terms vanish at r = 0. Hence, in Eq. (I22t AGx(O) cancels the 
n = term of AGx(r). The following series expansion is obtained for Cx(r), in view of Eqs. (|23, (|25j, and O: 

C (r) - V ^ ' 



(27r)'i/2 ^^ r(n+l)r(n+|) 



In Ught of Eq. i26\ . the sufficient condition i24\ is equivalent to the following 



kc 







f](-l)"+iii„(r) <c, (27) 



71=1 



where w„(r) are non-negative functions given by 



7-2" llogrP+'' 
ujr) = ——77^ :rMd) (28) 

and An{0) represents the following integral: 

A„(0)- J dke-+''+^f{k;e"), (29) 



and f{k; 0") is given by Eq. (|9}. The condition ( I27> is satisfied if the alternating series l)"+^u„(r) converges. 

An alternating series converges if it is absolutely convergent [39, p. 18]. According to the comparison test [39, 
p. 20], the series is absolutely convergent if m„ < Gm„, where ^ u„ is a convergent series, and G is a constant 
independent of n. 



Febraaiy 1, 2008 



DRAFT 



SUBMITTED TO IEEE TRANS. INFORM. THEORY 



11 



Using the mean value theorem [39, p. 65], the integral An{9) is evaluated as follows: 



Anie) f{kn;e") Jdkk 



2n+d+l 





/ ]^ 2n+2+d 



where /c„ g [0, kc], Vn. Let us define as f{k*; 6") — lim„_+oo/(fcn; d") the upper limit of the sequence /(fc„; 0"). 
The upper limit exists and is a finite number, since Vn, f{kn\9") < max{/(fc;0"), k E [0,kc]}. Then, using 
a = 2n + 2 + d, the following inequahty is obtained, V C G M : C > /(fc*; 6") 

A^{e)<f{k*;e")(^)<c(^). (31) 



Based on the inequality J31> . the sequence of absolute values, Un{r), of the initial series is bounded by the 
sequence Cu„, where: 

- ^ Voi^cO'' kc^|logr|i+^(kcr)2" 

{2Tr)dn (2n + 2 + d)r(n + l)r(n+|)' ^ ' 

To determine the convergence of the series X^J^Li we use d' Alembert's ratio test [39, p. 22], which states that 
the series converges absolutely if there is a fixed no, such that for all n > uq, |{t„+i/u„| < cq , where < cq < 1. 
Based on Eq. (I32> . the respective ratio is given by 

= (ker)2/3(Ti,d), (33) 



where 

(n + l + |)r(n + l)r(n+|) 



(3in,d) : 

^ (34) 



(n+l)(n + l + f) 

The function (3{n, d) is monotonically decreasing with n. For fixed kc, r let us define as no the smallest integer for 
which P{n, d) < (kcr)^-^. Then, \un+i/un\ < 1, Vn > no. This concludes the proof of sample path differentiability. 

■ 

Remark 4: Note that higher values of kc lead to higher threshold values no, implying a slower convergence of 
the series ( l26l and thus rougher SSRFs. Hence, kc provides a handle that permits controlling the roughness of the 
SSRF. 
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V. CovARiANCE OF One-dimensional FGC-SSRF Model 

The ID SSRF model can be applied to the analysis of time series and spatial data from one-dimensional samples 
(e.g., from drilling wells). 

Based on Eq. (|8|l, the ID covariance spectral density is given by the following expression 



G^,,{k;e)- ^^^^^ (fcO^ + (fcC)^- 
The covariance function is then obtained from the inverse Fourier transform of Eq. (|6j, i.e., 



Gx(r) 



dk ~ 

— G^{k) exp{jkr) 
Zn 



(35) 



= /dfc 



|QA(fc)| cos(fcr) 

T+mMFT{kW 

u 

Using the change of variables x = fc ^, and focusing on the boxcar kernel, we obtain 



Gx(r) 



Vo 



cos(x^ ^r) 



(36) 



1 + 771 x2 + x4 ■ 



Next, we calculate the variance and the integral scale of the covariance function for general kc, and we provide 
explicit asymptotic expressions for the covariance function for kc^ — > cx). It can be shown numerically that the 
asymptotic expressions are accurate for kc^ > 2, except for the differentiability at the origin. 
First, we define the dimensionless constants 



|2Tm 



1/2 



(37) 



Wl,2 



hiTAI 



1/2 



A. The Variance 

The variance is calculated based on Eq. ( I12> . 

Proposition 1 (The FGC-SSRF Variance): The variance is linearly proportional to 770, i.e.. 



where the function Vi (771 , x) is given by the following expressions, depending on the value of 771 : 



r 1 , /x2 + 2/3ix+l 
In 



4 f3i "' - 2 /3ix + 1 
'K 2^ 



Vi(77i,x) = < 



2/3: 



tan 



P2 



tan ^ (x) + — 



E 



(-1) 



i+i 



1=1,2 



tan 



for Iryil < 2, 
for 771 = 2 
for 771 > 2. 



(38) 



(39) 



(40) 
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Proof: The proof is given in Appendix I. ■ 

The variance of the Spartan model is a function of three parameters, namely 770, ?7i and kc^. Figure ^ displays 
the dependence of Vi(?7i,kc^) on kc^ for different values of the shape parameter rji in the range between —1.9 
and 4. For fixed 771, Vi(r7i, kc^) approaches the respective asymptotic limit for kci^ > 2. For fixed kc^, the function 
yii^lii^cC)! ^iid consequently the variance, decrease monotonically with increasing 771. 





vi 

---V2 - 









l# ' ' ' ' ' 1 

1 2 3 4 5 6 

X 



Fig. 1. Dependence of the function Vi(»)i,x) on x = kc^ for five different values of Tyi. 



B. The Integral Scale 

According to Eq. (I14> . the integral scale in d = 1 is given by 

770 27r^ 



(41) 



■G,(0) Vi(77i,kee)' 

A distinct feature of the SSRF covariance functions is the nonlinear dependence of the integral scale on the 
characteristic length ^ for kc^ < 2. However, if kc^ > 2, the integral scale Ii{9') becomes practically independent 
of the cutoff. Then, Ii{9') is essentially a function of only two variables: the shape parameter rji and the length ^. 
The dependence on <^ is linear in this asymptotic regime. More precisely, 

'4^/32, for|77i|<2, 

h{m,0=U^, for 7^1 = 2 (42) 

,2^ (wi + CJ2), for rji > 2. 

While the variance and the integral scale tend to asymptotic values for kci^ > 2, this behavior does not extend 
to the SSRF derivatives, i.e., to the integrals in Eq. ( I18> . which fail to converge with increasing kc. 
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C. Infinite-Band Covariance 

The covariance function can be evaluated explicitly for any combination of model parameters by means of the 
hyperbolic sine and cosine functions, as well as the sine and the cosine integrals. However, the resulting expressions 
are quite lengthy. Shorter asymptotic expressions are obtained, which are valid for kc^ > 2. More specifically: 



Proposition 2 (FGC-SSRF Covariance): The Spartan covariance depends linearly on the scale factor rjQ. For 
> 2 it becomes a function of the normalized distance h = |r|/^ and rji as follows: 



where the function Wi{h, rji) is given by the following: 



(43) 



Wi 



cos(/i/3i) ^ sin(/i/3i) 



4/32 



4/3i 



_\{l + h) 



1 /e^'^'^i 



Proof: The proof is given in the Appendix II. 



, for|r,i|<2, 
for 771 = 2 
for 771 > 2. 



(44) 



Corollary 2 (The auto-correlation function): The auto-correlation function is given by the equation: 



p(r) 



(^2 

cos(hPi) + — sin(/i/3i) 
Pi 



(l + /i)e-'\ 



for I771I <2, 
for r]i — 2 
for ?7i > 2. 



(45) 



LU2 - OJl 

Proof: By definition, the autocorrelation function is given by Px{r) Gx{r)/<j'^. The Eq. i45i follows from 
Eqs. (|39j, gOl, and (|44|l. ■ 

The autocorrelation function for 771 = 2 corresponds to the Whittle-Mattern function, Pi,{r) = r'^ Ki,{r) 
with = 3/2 [33]. It is interesting to note that the Whittle-Mattern covariance functions are obtained by solving 
a stochastic equation with an external white-noise forcing, while the Spartan covariance function is obtained from 
a Gibbs energy functional. For 7/1 > 2, an empirical correlation function is obtained [7]. The correlation function 
for 1 771 1 < 2 provides a class of positive definite functions in d = 1, which, to our knowledge, is new. 

Remark 5: The correlation function obtained for I771I < 2 in Eq. (145 > is not merely a superposition of permissible 
models, since the second term contains a sine function. However, the superposition of the two terms with the precise 
coefficients ensures the permissibility and the differentiability of the correlation function, in spite of the exponential 
term. 

Plots of the correlation function are shown in Figure (|3 for different values of the shape parameter For r/i < 
the correlation function oscillates, and the number of oscillations increases as 771 —2. The oscillations disappear 
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for positive values of ryi, and a monotonic decline of the correlations due to the exponential terms sets in. 




Fig. 2. Dependence of the autocorrelation function on distance for different values of the shape parameter r]i with fixed ^ = 0.1. 



VI. COVARIANCE OF THREE-DIMENSIONAL FGC-SSRF MODEL 

The spectral representation, i.e., Eq. ilQi . of the isotropic Spartan co variance in d = 3 is given by 

Voe J^,^ k^/'Ji/2ikr)\Qx{k)\' 
'^^'■^ " (2^)3/2 ri/2 J ""1 + rj,{k^)^ + 



Using the identity 



it follows that 



/2y/^sin(r) 



VaS,^ f „ k sm{kr) \ Qx{k) 



|2 



i + m{kO^ + {kiy 



Using the transformation u ~ k£^, and the boxcar kernel spectral density, we find 

n I \ ^?o^ /" J u sin(r<~i) 


A. The Variance 

The SSRF variance is calculated based on Eq. M2\ . We use the dimensionless quantities /3i, /32, wi, and ijJ2 
defined in Section dVt . 
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Proposition 3 (FGC-SSRF Variance): The variance of the Spartan co variance in d = 3 is given by 



Gx(0) = ^V3(r,i,kee) 

where 



(48) 



f 1 , /^x2-2/3ix+l 
In 



4/3i \x^ + 2f3ix+l 



V3 = < 



2/32 ^ 



tan 



1^ 



tan~^(x) 
2 



1+X2' 



1 = 1,2 

Proof: The proof is presented in the Appendix III. 



for bnl <2, 

for 771 = 2 
for ?7i > 2. 



(49) 




Fig. 3. Dependence of the function V3(r;i,x) on x = kc^ for five different values of 771. 



B. The Integral Scale 

According to Eq. (I14t . the integral scale in d = 3 is given by 

As seen in Figure Q, the approach of V3(77i,kc^) to the asymptotic limit is slower than in d = 1. The integral 
scale 13(6') becomes practically independent of the cutoff for kc^ > 5. The asymptotic expressions of the integral 



Vo 



Gx(0)_ 



(50) 
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scale are: 



Tt{ll!i + LU2) 



2C 



1/3 



for I771I < 2, 
for 771 = 2 
for ?7i > 2. 



(51) 



C. Infinite-Band Covariance 

As in d = 1, the Spartan covariance function can be evaluated explicitly for any 6 by means of the hyperbolic sine 
and cosine functions, as well as the sine and the cosine integrals. Here we give the asymptotic (in kc) expressions: 



Proposition 4 (FGC-SSRF Covariance): The covariance for kc ^ 00 is expressed as follows: 

Gx(r) = ^M/3(r/e,r7i), 

ZTT 

where the function W3(/i, 771), h = r/^ is given by the following: 



-ft/3i 



sin ( ft/32 )' 



1 



2A 



for bnl <2, 
for rji ~ 2 
for 771 > 2. 



(52) 



(53) 



Proof: The proof is given in the Appendix IV. ■ 
The known exponential covariance is obtained for rji = 2, while for |?7i| < 2 a product of two permissible 
models, i.e., exp(/i) and sm{h)/h, is obtained. The covariance model obtained for 771 > 2 is new, at least to our 
knowledge. 



Corollary 3 (The autocorrelation fimction): The auto-correlation function is given by the equation: 

sin (/1/32) 



p(r) 



hP2 



, for |?7i| < 2, 
for 7/1 = 2 
for 771 > 2. 



(54) 



. h{uj2 — LUl) 

Proof: It follows from the definition of the autocorrelation function, as well as Eqs. ( 1?^ . and ( 15^ . ■ 

The dependence of the autocorrelation function on distance for various values of 771 is shown in Figure 0. Note 
that the negative hole for r/i = — 1 is significantly less pronounced compared to the d — 1 case. 

Remark 6: The Spartan covariances obtained for infinite band in d = 3 are continuous but non-differentiable, in 
contrast with the d = 1 case. 
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Fig. 4. Dependence of the autocorrelation function on distance for different values of the shape parameter rji with fixed ^ = 0.1 



In d = 2 the integral of the covariance function can not be evaluated explicitly, even in the infinite band case. 
However, the covariance functions obtained for d = 3 are also permissible d = 1,2 [1, pp. 30-31]. Of course, they 
are not derived from the FGC-SSRF model in d < 3. 

VII. Discussion and Conclusions 

We showed that FGC Spartan Spatial Random Fields are analytic in any dimension, provided that the covariance 
spectral density is band-limited. We calculated the variance and the integral scale of the FGC-SSRF covariance 
functions in one and three dimensions. We also obtained explicit expressions for the covariance function at the 
asymptotic limit (i.e., the infinite-band limit). Depending on the value of the shape parameter 771, the resulting 
expressions were shown to recover known models or to yield new covariance functions. 

Explicit expressions for the covariance functions in d = 1,3 are also possible in the pre-asymptotic limit. 
However, they are not given here since they are very lengthy, and in practice it may be preferable to integrate 
Eq. ilO\ numerically. It was also shown that the asymptotic expressions are quite accurate for the covariance 
function (but not for its derivatives) when the product kc^ exceeds a finite, dimensionality-dependent threshold. 
This result has practical applications in SSRF parameter inference, since the procedure used involves matching 
ensemble constraints, which are expressed in terms of the covariance function, with respective sample constraints 
[18], [35]. 

Explicit expressions for the covariance function were not found in d = 2, where the presence of Jo{kr) prohibits 
closed form integration in Eq. dlOt . Explicit expressions for the variance are given in [22], and for the integral scale 
in [35]. 
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The FGC-SSRF covariance functions permit a continuous transition between smooth (analytic) and rough (non- 
analytic) states, by controlling the spectral band cutoff kc. This property is shared by the Whittle-Matern covariance 
functions. 

Besides providing new covariance functions, the SSRF idea focuses on representing spatial structure using energy 
functionals with clear physical interpretation. These lead to new possibilities for parameter inference and spatial 
interpolation in geostatistical appUcations, which are the target of continuing investigations [18], [23], [35], [37], 

Appendix I 
Proof of Proposition 1 

Proof: We calculate the integrals required for evaluating the variance of the ID SSRF. 



Step (i): |77i| < 2. We define the characteristic polynomial n(x) ^ x"' + ?7ix2 ^ ^ expand n(a;) as follows: 

n(x) = (x^ + 2/3ix + l)(x2 - 2/3ix + 1), 

where /3i is given by Eq. (I37> . 

Using partial fraction expansion we obtain 

1 1 / 2x + 2/3i 



2x - 2/3i 

n(x) ' 8/3i + 2/3ix + 1 " x2 - 2/3ix + 1 



+ - 



1 



1 



1 



2 + ?7i 

where /32 is given by Eq. (I38> . Then, by direct integration the following is obtained: 

1 



^_ _ J_ 1^ x2 +2/3ix+ 1 
n(x) " 8/3i ^^U'-2/3ix+l 



1 

W2 



tan 



x + A 
P2 



tan 



1 ( X - /3i 

/32 



(55) 



Step (ii): 771 = 2. In this case one obtains: 

1 1 l+x2 1 l-x2 



n(x) 2 (1+X2)2 2 (1+X2)2 

1 l-x2 



2(1 + x2) 2(1 +x2) 

It follows that 

1 tan~^(x) X 



(56) 



c?x ■ 



n(x) 2 2(1 + x2)' 

Step (iii): 771 > 2. In this case n(x) is expressed as follows: 

n(x) = {yi^ +ujI){i^ 
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where uJi,2 are given by Eq. ( I38> . We find 

1 



dx 



n(x) 



-Idx 



1 



1 



1 

A 



1 

UJi 



-tan 



1 



-tan 



LU2 



X 

UJ2 



Appendix II 
Proof of Proposition 2 

Proof: Let us define the normalized distance h = |r|/^. We evaluate the integrals using integration in the 
complex plane. More specifically, based on the residue theorem and Jordan's lemma [8, pp. 358-370], we obtain 



oo 

I 



dx 



cos(/ix) 



= i TTj ^ Res^ 



exp(jhz) 



(57) 



where ^Res"*" denotes the sum of the residues of the function Q{z) in the upper half plane, and 9?[A(z)] 

denotes the real part of the complex function A{z). 



Step (i): |?7i| < 2. 

Since n(z) — (z^ + 2Piz + l) (z^ — 2(3iz + l) the imaginary poles in the upper half plane are z± — ±/3i +j /?2- 
Thus, we obtain 



exp(j hz) 

n(z) 



Z — Z± 



J 

2/32 



and 



77i-2 + jA 7?i-2-jA 
exp(j hz) 



n(z) 



-/1/32 



Z — Z± 

(771 - 2) cos (hPi) - A sin {hPi) 



(4771 - 8)/32 



Step (ii): r/i = 2. 

Since n(z) = (1 + z^)^, in the upper-half plane there is a double imaginary pole, z+ ~ j. Hence, 

oc 


7r(l + fe) 
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Step (iii): rji > 2. Using the variables u)j, recall from Appendix I that 



1 1 



Thus 



oo 

/ 



rfx 



1 



1 



cos(x/i) 



Appendix III 
Proof of Proposition 3 

Proof: We calculate the integrals required for evaluating the variance of the 3D SSRF. 



Step (i): < 2. 
Using the partial fraction expansion we obtain 

X2 1 / 2X - 2; ii 



2x + 2/^1 



n(x) 8/3i V x2 - 2/3ix + 1 x2 + 2/3ix + 1 

1 x^ + l 

^2 (x2 + 2/3ix+l)(x2-2/3ix+l)" 



By direct integration we obtain 



^^^nix) "8/3/°^U' + 2/3ix+l 



+ 7^3" — 3 — + tan — - — 
4/32 I \ (^2 J V /32 



Step (ii): r?i = 2. In this case the partial fraction expansion becomes 



It follows that 



Step (iii): 771 > 2. 
In this case 



/ 



dx 



n(x) 



1 



l-x^ 



n(x) 2(1 +X2) 2(l+x2)2 



tan-i(x)-- 



+ X2 



n(x) = (x2 + UjI) (x2 + UjI 



February 1, 2008 



DRAFT 



SUBMITTED TO ffiEE TRANS. INFORM. THEORY 



22 



and 



It follows that 



n(x) AVx2+w| x2+w2y 



/ 



dx 



n(x) 



1 

A 



UJ2 tan ^ I — ) — wi tan ^ ( — 
W2/ V'^i 



Appendix IV 
Proof of Proposition 4 

Proof: We follow the same procedure for integrating the covariance spectral density as in the Appendix II. 
In = 3 the respective integral is given by: 



/ 



n(x) 



TTj^Res" 



-3 hz 



Uiz) 



where X^Res denotes the sum of the residues of the function Q{z) in the lower half plane, and 

denotes the imaginary part of the complex function A{z). 

Step (i): |77i| < 2. 
The imaginary poles in the lower half plane are 



We have 



Thus, 



z e 



-]hz 



U{z) 



J 
2/32 



Z — Z± 

2-r?i 2-r/i . 

= ^ sin(/i/3i) e-'''^\ 



OO 

/ 



rfx 



(58) 



Step (ii): 771 = 2. 



(ix 



x sin(x h) 

n(x) 



00 

/X sin(x h) 7T 
'^'^ (X2 + 1)2 =- 
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Step (iii): 771 > 2. Using n(z) — (z^ + wf) (z^ + luI), we obtain 

zsin(z/i) 1 (zsm(zh) zsm(zh) 



n(z) 



The poles in the lower half plane are z+ — —j L02 and z_ 



xsin(x/i) _ ^ I 

'^'^ n(x) - A ^ 1 



z^ 

J wi . It then follows that; 



j Res 



"z e"^'*^' 




_z2 + a;|_ 


•J 







TT 

2A 



Z^+Lof 
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